1 Introduction

In this notebook we will perform the gene-peak links filtering by deferentially accessible regions (DARs) and, then, by the differentially expressed genes (DEGs). We need to remember that we are only analyzing BCL6 and PRDM1 genes links as proof of concepts and to reduce the computational time. Also, you will see how the number or link significantly decrease and how they look at the coverage plot.

In this report we will continue working with the previous Seurat object saved as 4.tonsil_bcell_linkpeaks.rdsin the results/R_object directory.

Load packages 
BiocManager::install("Repitools")
knitr::opts_chunk$set(echo = TRUE)
library(Signac)
library(Seurat)
library(ggplot2)
library(ggpubr)
library(tidyverse)
library(EnsDb.Hsapiens.v86)
library(plyr)
library(reshape2)
library(data.table)
library(GenomicRanges)
library(harmony)
library(hdf5r)
library(stringr)
library(ggpubr)
library(RColorBrewer)
library(magick)
library(knitr) 
library(biovizBase)
library(patchwork)
library(Repitools)
library(qlcMatrix)
library(kableExtra)
library(BSgenome.Hsapiens.UCSC.hg38)
set.seed(123)

2 Loading data

2.1 Directory path

# Paths
path_to_objects <- here::here("results/R_objects/")
path_to_tables <- here::here("results/tables/")

2.2 Load data

tonsil_bcell <- readRDS(paste0(path_to_objects,"4.tonsil_bcell_linkpeaks.rds"))
deg<- read_csv(paste0(path_to_tables,"tonsil_bcell_deg.csv"))
## New names:
## Rows: 3696 Columns: 8
## ── Column specification
## ──────────────────────────────────────────────────────── Delimiter: "," chr
## (3): ...1, cluster, gene dbl (5): p_val, avg_log2FC, pct.1, pct.2, p_val_adj
## ℹ Use `spec()` to retrieve the full column specification for this data. ℹ
## Specify the column types or set `show_col_types = FALSE` to quiet this message.
## • `` -> `...1`
dars<- read_csv(paste0(path_to_tables,"tonsil_bcell_dars.csv"))
## New names:
## Rows: 25861 Columns: 8
## ── Column specification
## ──────────────────────────────────────────────────────── Delimiter: "," chr
## (3): ...1, cluster, gene dbl (5): p_val, avg_log2FC, pct.1, pct.2, p_val_adj
## ℹ Use `spec()` to retrieve the full column specification for this data. ℹ
## Specify the column types or set `show_col_types = FALSE` to quiet this message.
## • `` -> `...1`

3 Joint UMAP

Idents(tonsil_bcell)<-"annotation_level_1"
DimPlot(tonsil_bcell,label = TRUE, reduction = "wnn.umap", pt.size = 0.5,cols = c("#a6cee3", "#1f78b4","#b2df8a", "#fb9a99","#e31a1c")) + ggtitle("Joint UMAP")+NoLegend()

5 Coverage plot

coverage_extend{r}is a function that create a coverageplot for a certain gene and region using different upstream and downstream extension region. 

coverage_extend <- function(x,y,seuratobject){purrr::map(y, function(y) {
  
  p <- CoveragePlot(
  object = seuratobject,
  region = x,
  features = x,
  expression.assay = "RNA",
  idents = idents.plot,
  extend.upstream = y,
  extend.downstream = y
  #tile = TRUE
  )
  p & scale_fill_manual(values = cols_cluster)
})}
ranges.show <- StringToGRanges("chr3-187721377-187745725")

# set the colors will have  each chromatin accesibility  profile of each cell type. That match with the colors in the UMAP
cols_cluster <- c("#a6cee3", "#1f78b4","#b2df8a", "#fb9a99")
idents.plot <- c("GCBC", "NBC", "MBC", "PC")

8 Save

Session Info 
sessionInfo()
## R version 4.2.1 (2022-06-23)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur ... 10.16
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRlapack.dylib
## 
## locale:
## [1] es_ES.UTF-8/es_ES.UTF-8/es_ES.UTF-8/C/es_ES.UTF-8/es_ES.UTF-8
## 
## attached base packages:
## [1] stats4    stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] BSgenome.Hsapiens.UCSC.hg38_1.4.4 BSgenome_1.64.0                  
##  [3] rtracklayer_1.56.1                Biostrings_2.66.0                
##  [5] XVector_0.38.0                    kableExtra_1.3.4                 
##  [7] qlcMatrix_0.9.7                   sparsesvd_0.2-1                  
##  [9] slam_0.1-50                       Matrix_1.5-3                     
## [11] Repitools_1.42.0                  patchwork_1.1.2                  
## [13] biovizBase_1.44.0                 knitr_1.40                       
## [15] magick_2.7.3                      RColorBrewer_1.1-3               
## [17] hdf5r_1.3.7                       harmony_0.1.1                    
## [19] Rcpp_1.0.9                        data.table_1.14.4                
## [21] reshape2_1.4.4                    plyr_1.8.8                       
## [23] EnsDb.Hsapiens.v86_2.99.0         ensembldb_2.20.2                 
## [25] AnnotationFilter_1.20.0           GenomicFeatures_1.48.4           
## [27] AnnotationDbi_1.60.0              Biobase_2.58.0                   
## [29] GenomicRanges_1.50.1              GenomeInfoDb_1.34.2              
## [31] IRanges_2.32.0                    S4Vectors_0.36.0                 
## [33] BiocGenerics_0.44.0               forcats_0.5.2                    
## [35] stringr_1.4.1                     dplyr_1.0.10                     
## [37] purrr_0.3.5                       readr_2.1.3                      
## [39] tidyr_1.2.1                       tibble_3.1.8                     
## [41] tidyverse_1.3.2                   ggpubr_0.5.0                     
## [43] ggplot2_3.4.0                     SeuratObject_4.1.3               
## [45] Seurat_4.2.1                      Signac_1.8.0                     
## [47] BiocStyle_2.24.0                 
## 
## loaded via a namespace (and not attached):
##   [1] Hmisc_4.7-1                 ica_1.0-3                  
##   [3] svglite_2.1.0               RcppRoll_0.3.0             
##   [5] Rsamtools_2.12.0            lmtest_0.9-40              
##   [7] rprojroot_2.0.3             crayon_1.5.2               
##   [9] MASS_7.3-58.1               nlme_3.1-160               
##  [11] backports_1.4.1             reprex_2.0.2               
##  [13] rlang_1.0.6                 ROCR_1.0-11                
##  [15] readxl_1.4.1                irlba_2.3.5.1              
##  [17] limma_3.52.4                filelock_1.0.2             
##  [19] BiocParallel_1.30.4         rjson_0.2.21               
##  [21] bit64_4.0.5                 glue_1.6.2                 
##  [23] sctransform_0.3.5           parallel_4.2.1             
##  [25] spatstat.sparse_3.0-0       vsn_3.64.0                 
##  [27] spatstat.geom_3.0-3         haven_2.5.1                
##  [29] tidyselect_1.2.0            SummarizedExperiment_1.28.0
##  [31] fitdistrplus_1.1-8          XML_3.99-0.12              
##  [33] zoo_1.8-11                  GenomicAlignments_1.32.1   
##  [35] xtable_1.8-4                magrittr_2.0.3             
##  [37] evaluate_0.18               cli_3.4.1                  
##  [39] zlibbioc_1.44.0             rstudioapi_0.14            
##  [41] miniUI_0.1.1.1              sp_1.5-1                   
##  [43] bslib_0.4.1                 rpart_4.1.19               
##  [45] fastmatch_1.1-3             shiny_1.7.3                
##  [47] xfun_0.35                   cluster_2.1.4              
##  [49] caTools_1.18.2              KEGGREST_1.38.0            
##  [51] ggrepel_0.9.2               listenv_0.8.0              
##  [53] png_0.1-7                   future_1.29.0              
##  [55] withr_2.5.0                 ggforce_0.4.1              
##  [57] bitops_1.0-7                cellranger_1.1.0           
##  [59] pillar_1.8.1                gplots_3.1.3               
##  [61] cachem_1.0.6                fs_1.5.2                   
##  [63] vctrs_0.5.1                 ellipsis_0.3.2             
##  [65] generics_0.1.3              tools_4.2.1                
##  [67] foreign_0.8-83              tweenr_2.0.2               
##  [69] munsell_0.5.0               DelayedArray_0.24.0        
##  [71] fastmap_1.1.0               compiler_4.2.1             
##  [73] abind_1.4-5                 httpuv_1.6.6               
##  [75] Ringo_1.60.0                plotly_4.10.1              
##  [77] rgeos_0.5-9                 GenomeInfoDbData_1.2.9     
##  [79] gridExtra_2.3               DNAcopy_1.70.0             
##  [81] edgeR_3.38.4                lattice_0.20-45            
##  [83] deldir_1.0-6                utf8_1.2.2                 
##  [85] later_1.3.0                 BiocFileCache_2.6.0        
##  [87] jsonlite_1.8.3              affy_1.74.0                
##  [89] scales_1.2.1                docopt_0.7.1               
##  [91] pbapply_1.6-0               carData_3.0-5              
##  [93] genefilter_1.78.0           lazyeval_0.2.2             
##  [95] promises_1.2.0.1            car_3.1-1                  
##  [97] latticeExtra_0.6-30         goftest_1.2-3              
##  [99] spatstat.utils_3.0-1        reticulate_1.26            
## [101] checkmate_2.1.0             rmarkdown_2.18             
## [103] cowplot_1.1.1               webshot_0.5.4              
## [105] Rtsne_0.16                  dichromat_2.0-0.1          
## [107] uwot_0.1.14                 igraph_1.3.5               
## [109] survival_3.4-0              yaml_2.3.6                 
## [111] systemfonts_1.0.4           htmltools_0.5.3            
## [113] memoise_2.0.1               VariantAnnotation_1.42.1   
## [115] BiocIO_1.6.0                locfit_1.5-9.6             
## [117] here_1.0.1                  gsmoothr_0.1.7             
## [119] viridisLite_0.4.1           digest_0.6.30              
## [121] assertthat_0.2.1            mime_0.12                  
## [123] rappdirs_0.3.3              RSQLite_2.2.18             
## [125] future.apply_1.10.0         blob_1.2.3                 
## [127] preprocessCore_1.58.0       splines_4.2.1              
## [129] Formula_1.2-4               labeling_0.4.2             
## [131] googledrive_2.0.0           ProtGenerics_1.28.0        
## [133] RCurl_1.98-1.9              broom_1.0.1                
## [135] hms_1.1.2                   modelr_0.1.10              
## [137] colorspace_2.0-3            base64enc_0.1-3            
## [139] BiocManager_1.30.19         nnet_7.3-18                
## [141] sass_0.4.2                  bookdown_0.30              
## [143] RANN_2.6.1                  fansi_1.0.3                
## [145] tzdb_0.3.0                  truncnorm_1.0-8            
## [147] parallelly_1.32.1           R6_2.5.1                   
## [149] grid_4.2.1                  ggridges_0.5.4             
## [151] lifecycle_1.0.3             curl_4.3.3                 
## [153] ggsignif_0.6.4              googlesheets4_1.0.1        
## [155] affyio_1.66.0               leiden_0.4.3               
## [157] jquerylib_0.1.4             RcppAnnoy_0.0.20           
## [159] spatstat.explore_3.0-5      htmlwidgets_1.5.4          
## [161] polyclip_1.10-4             biomaRt_2.52.0             
## [163] crosstalk_1.2.0             timechange_0.1.1           
## [165] rvest_1.0.3                 globals_0.16.1             
## [167] htmlTable_2.4.1             spatstat.random_3.0-1      
## [169] progressr_0.11.0            codetools_0.2-18           
## [171] matrixStats_0.62.0          lubridate_1.9.0            
## [173] gtools_3.9.3                prettyunits_1.1.1          
## [175] dbplyr_2.2.1                gtable_0.3.1               
## [177] DBI_1.1.3                   highr_0.9                  
## [179] tensor_1.5                  httr_1.4.4                 
## [181] KernSmooth_2.23-20          stringi_1.7.8              
## [183] vroom_1.6.0                 progress_1.2.2             
## [185] farver_2.1.1                annotate_1.74.0            
## [187] DT_0.26                     xml2_1.3.3                 
## [189] restfulr_0.0.15             interp_1.1-3               
## [191] scattermore_0.8             bit_4.0.5                  
## [193] jpeg_0.1-9                  MatrixGenerics_1.10.0      
## [195] spatstat.data_3.0-0         pkgconfig_2.0.3            
## [197] gargle_1.2.1                rstatix_0.7.1              
## [199] Rsolnp_1.16